arXiv:1509.07353vl [cond-matsoft] 24 Sep 2015 


Crack patterns over uneven substrates 


Pawan Nandakishore and Lucas Goehrin^ 

Max Planck Institute for Dynamics and Self-Organization (MPIDS), 37077 Gottingen, Germany 

(Dated: September 25, 2015) 

Cracks in thin layers are influenced by what lies beneath them. From buried craters to crocodile 
skin, crack patterns are found over an enormous range of length scales. Regardless of absolute size, 
their substrates can dramatically influence how cracks form, guiding them in some cases, or shielding 
regions from them in others. Here we investigate how a substrate’s shape affects the appearance of 
cracks above it, by preparing mud cracks over sinusoidally varying surfaces. We find that as the 
thickness of the cracking layer increases, the observed crack patterns change from wavy to ladder¬ 
like to isotropic. Two order parameters are introduced to measure the relative alignment of these 
crack networks, and, along with Fourier methods, are used to characterise the transitions between 
crack pattern types. Finally, we explain these results with a model, based on the Griffith criteria of 
fracture, that identifies the conditions for which straight or wavy cracks will be seen, and predicts 
how well-ordered the cracks will be. Our metrics and results can be applied to any situation where 
connected networks of cracks are expected, or found. 


I. INTRODUCTION 


From networks of cracks across the surfaces of planets 
to artistic craquelure in pottery glazes and paint¬ 
ings m\, or bespoke fracture textures in thin films El El, 
crack patterns naturally occur on a wide variety of scales. 
They are even implicated in the formation of scales on 
crocodiles m- All these structures represent the out¬ 
come of failure processes that occur in a brittle medium 
supported by an uneven substrate. Here we show how 
the shape of the substrate plays an important role in de¬ 
termining the geometry of the resulting cracks. We will 
demonstrate how to use the substrate to template frac¬ 
tures, explore the conditions necessary to generate a vari¬ 
ety of patterns, and introduce general metrics to quantify 
the orientational ordering of any set of connected cracks. 

To introduce the concerns of our work we will first 
briefly present three examples of cracks on the scales of 
kilometres, centimetres, and micrometers. These motiva¬ 
tional examples, and those shown additionally in Fig. 
also aim to highlight the need for quantitative approaches 
to the analysis of fracture patterns. 

On large length scales, radial graben are associated 
with buried craters on Mercury [ami and Mars BE]- 
These features are due to the repeated filling and cooling 
of lava in impact craters. The uneven stresses caused 
by the edge of the crater are believed to lead to wrinkle 
ridges above the crater rim [4]. Contraction cracks also 
form within the interior of the crater, leading to a mixture 
of polygonal and circular graben patterns [3[4]. The 
radial symmetry of the substrate imposes a characteristic 
structure to the such cracks, as in Fig. [^a). 

On a more intermediate length scale, craquelure 
(Fig-Hb, c)) is the result of differential strains between a 
coating and its substrate iSHa E]. This may be due 
to tension in the framing of a painted canvas [5], for 
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FIG. 1: Examples of crack networks, oriented by different ef¬ 
fects. (a) Gracks in lava over a buried crater in Goethe Basin, 
Mercury HU. (b) Graquelure pattern of traditional oil-on- 
wood technique, (c) Paint cracks on a beam, reflecting the 
growth rings of the underlying wood, (d) Gracks in dried 
blood droplets, in a range of humidities, are affected by the 
way the droplet dries (courtesy of D. Brutin, adapted from 
HU, Fig. 6). (e) Skin cracks on embryonic crocodiles are the 
result of differential growth of the skin layers, and form tiles 
that mature into scales (courtesy of M. Milinkovitch, adapted 
from HU, Fig. lb), (f) Desiccation cracks in pastes are influ¬ 
enced by any vibration or flow prior to drying (courtesy of A. 
Nakahara, adapted from [13], Fig. 3b). 


example, or changes in the relative humidity at which 
a painting has been stored HU. Alternatively, it may 
result from an artist’s choice of a pottery glaze with a 
higher coefficient of thermal contraction than the clay 
it covers. Such craquelure patterns are extremely hard 
to forge, and can thus aid in the attribution of paint¬ 
ings |6| . They have been shown to be influenced by local 
painting style, choice of substrate, including, for paint on 
wood, even the species of wood used [6l [TSK The analy¬ 
sis of similar patterns, in dried blood (Fig. Wd)), is also 
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being developed as a forensic tool, for example to help 
determine the relative humidity when the blood was spilt 

[nms]. 

On the micrometer scale, recent attention has focussed 
on controlling cracks. For example, Nam et al. [8] 
showed how to selectively initiate, refract, and stop both 
straight and wavy cracks in silicon nitride thin films, 
vapour-deposited on silicon wafers, by adding notches, 
cutting stairs, or introducing interfacial layers between 
the wafer and film. Kim et al. [9] instead looked at soft 
substrates, with crack formation guided by strategically 
placed notches along grooves; their work aims to control 
failure in situations, such as flexible electronics, where 
it must be expected. They find the most regular results 
when they design for crack spacings that are near the 
natural crack spacing of the brittle film. 

A common feature of all these situations is that the 
cracking layer is thin, compared to its lateral extent. 
Nonetheless, the layer thickness imposes the natural 
length-scale of such problems, and the spacing of the 
cracks in each case reflects this thickness [ilTTHIS]. Here 
we study the model problem of contraction cracks in 
‘thin’ elastic layers adhered to rigid substrates and show 
how the shape of the substrate also plays an important 
role in determining the geometry of the crack pattern. 

For this end we choose the system of desiccation cracks 
in drying clays as a simple, yet representative, experi¬ 
ment. As clays dry they tend to shrink, uniformly. So- 
called channel cracks open in the film, to relieve its grow¬ 
ing strain energy (see e.g. [2Q]422] b The appropriate 
physical model for stress in a wet clay, poroelasticity, has 
an exact mathematical analogy with the thermoelastic 
model of how most solids shrink as they cool [23f[25] . 
In biology, comparable conditions can also result from 
the differential growth of layered tissues [261128] . This is 
known to lead to the formation of scales on the heads 
of Nile crocodiles [To] (Fig. [^e)), for example, while a 
similar mechanism has been proposed to account for leaf 
venation patterns [2i|27|. 


II. MATERIALS AND METHODS 

Briefly, bentonite clay was dried over rigid patterned 
substrates, where it cracked; a schematic summary of the 
experiment is given in Fig. |^a). Compared to the wet 
clay, the substrate is undeformable and unbreakable. 

Each substrate consisted of a 20x20 cm^ plate. Five si¬ 
nusoidal plates were cut from acrylic blocks by computer- 
numerical-control (CNC) milling, with a resolution of 
200-400 fim. These plates were corrugated; their sur¬ 
face 2 ; varied sinusoidally along one direction, x, with a 
constant profile along the other, y, direction such that 
2 : = Asin(27rx/A). The wavelengths A and amplitudes 
A of the plates are given in Table |I| Two pairs of plates 
had the same dimensionless ratios a = A/A, in order 
to check the scale-independence of the fracture patterns. 
Two additional plates were prepared by 3D printing from 
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FIG. 2: Sketch of the drying experiment, where (a) a clay 
slurry of average initial thickness H dries over a sinusoidally 
varying substrate, of amplitude A and wavelength A. Both 
(b) linear and (c) radial wave plates were used as substrates. 


TABLE I: Amplitudes A and wavelengths A of the corrugated 
substrates (wave plates) used to template crack patterns. Lin¬ 
ear wave plates are simply numbered, while radial wave plates 
are indicated by the suffix r. 


Plate 

A (cm) 

A (cm) 

a = A/ A 

1 

0.25 

2 

0.125 

2 

0.25 

1 

0.25 

3 

0.5 

2 

0.25 

4 

0.25 

0.5 

0.5 

5 

0.5 

1 

0.5 

Ir 

0.25 

1 

0.25 

2r 

0.25 

0.5 

0.5 


acrylic photopolymer, with a vertical resolution of 30 
/im, and a contour resolution of 5 jam (4D Concepts). 
These plates had radial wave patterns, with surface pro¬ 
files (see Table that varied periodically in the radial 
direction, r, away from the centre of the plate, such that 
^ = A cos(27rr/A). Einally, a flat acrylic plate was used 
for control experiments. All plates were smooth to the 
touch. Eor each plate 10 cm high acrylic walls were glued 
to each edge, to make a waterproof container. 

Slurries of potassium bentonite (Acros Organics K-10 
montmorillonite) were prepared by mixing a set weight 
of bentonite powder with twice as much water, by mass. 
Slurries were then stirred until smooth and immediately 
poured evenly into a container, over a substrate. At this 
concentration the slurries behaved as fluids, rather than 
pastes. These conditions avoid the memory effect [291130], 
whereby the flow or vibration of a yield-stress paste can 
affect the way it cracks (see Eig. Bf))- We checked 
that this effect was absent by drying slurries over flat 
substrates, and finding isotropic crack patterns. 
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After pouring, the bentonite settled slightly, separat¬ 
ing into a dense clay covered by a thin water layer. We 
measured the thicknesses of these submerged clay lay¬ 
ers in flat-bottomed containers after 6 hours of settling, 
without significant evaporation, and found that the con¬ 
solidated wet clay consistently contained 0.49 grams of 
bentonite per cubic centimetre. For slurries poured over 
the sinusoidal substrates we therefore used the total mass 
of the bentonite as a proxy for the average layer thickness 
H of the wet clay layer. All experiments were performed 
in the range of H from 3 to 15 mm. 

Once a slurry had been poured, it was left to dry in 
a quiet, warm environment for 6-48 hours, depending on 
the amount of clay used. When drying is complete, the 
bentonite becomes lighter in colour; no further cracks ap¬ 
pear after this time. For the linear wave plates, drying 
was accelerated by halogen heat lamps, which raised the 
surface temperature of the clay to ^ 50° C. As elevated 
temperatures would melt the printed plates, drying over 
the radial plates was enhanced by raising the room tem¬ 
perature to about 30° C. During drying, the appearance 
and evolution of the cracks were imaged every 5 min¬ 
utes using a digital SLR camera (Nikon D5000 series), 
mounted above the drying container. By comparing the 
initial wet height H to the dry height of various layers, 
we estimate a total vertical strain of 0.27, during drying. 


A. Image processing 

Images were pre-processed using Imaged and Matlab, 
as shown in Fig. Each image was cropped to isolate 
the cracked clay, then converted to grayscale by averaging 
across all colour channels, and contrast-enhanced. Next, 
a band-pass filter was applied to remove high-frequency 
noise (3-pixel cutoff to reduce single-pixel noise) and low- 
frequency variations in lighting intensity (40-pixel cutoff, 
comparable to the width of the largest cracks). Finally, 
each image was thresholded to give a black-and-white 
representation of the crack network, as in Fig. |^b). 

To characterise the crack patterns we will focus on the 
shape of the crack network, and on the relative positions 
and orientations of features such as the vertices between 
cracks. For some measurements, we thus chose to elimi¬ 
nate unconnected cracks and spurs from the binary im¬ 
ages. Isolated cracks, and other noise, were removed by 
deleting objects below a threshold size. Every image was 
then thickened, by 3 pixels, to ensure that small breaks in 
the crack network, or en-passant cracks m, were filled 
in. The result was then eroded to leave a skeleton con¬ 
sisting of one-pixel-wide cracks that connected various 
vertices or nodes, and which delimited isolated patches, 
or peds. Cracks that branched off the main network, but 
stopped without connecting to anything, were then re¬ 
moved. Both the binary images and the skeletons were 
used as the basis for further analyses. 

Einally, we note here that the error bars given through¬ 
out this work are estimated by subdividing the processed 


(a) (b) 



FIG. 3: Image processing. Each image (a) is cropped, 
contrast-adjusted, filtered and thresholded to (b) yield a 
black-and-white description of the crack positions. The shape 
of the substrate (plate I) is sketched below the images. 



FIG. 4: Evolution of crack patterns: (a) for thin layers, wavy 
cracks are seen; (b) as the layer height increases straight 
cracks are observed, with laddering cracks between them; and 
(c) for thick layers the cracks are isotropic. The substrate 
(plate 3) and clay thicknesses are sketched below each panel. 


images into quarters, analysing each subsection sepa¬ 
rately, and then calculating the standard deviation of the 
four different results. 


III. RESULTS AND DISCUSSION 

We found that as the thickness of the brittle layer 
increases, its contraction cracks slowly lose information 
about their substrate. However, in order to understand 
the effect of an uneven substrate, we must first briefly 
summarise the base case of cracks in flat bentonite layers. 
Eor the flat control plate, isotropic crack patterns were 
seen for all layer thicknesses {H from 4-15 mm). Cracks 
appeared one by one, sequentially breaking the clay into 
smaller and smaller pieces, as in [32l|33]. This process 
stops at some finite spacing, when the pattern saturates 
[T8l[25] - further drying increases the widths of the cracks 
in these ‘mature’ networks, as in [saiss], rather than cre¬ 
ate new ones. We confirmed that the natural crack spac¬ 
ing over the flat plate is between 2-3 times the wet layer 
thickness iL, using the line-dropping method [36]. The 
results are compatible with the crack spacings reported 
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FIG. 5: A phase diagram showing where wavy, ladder-like 
and isotropic crack patterns were seen, for the linear wave 
plates. 


in [36], but consistently smaller, as the dry heights of 
bentonite layers were reported in that study. 

For clays dried over the sinusoidally varying substrates 
the typical sequence of crack patterns, after drying is 
complete, is shown in Fig. A wavy pattern developed 
in thin layers, with the first cracks meandering back and 
forth as they travelled along the troughs of the underly¬ 
ing substrate. Various shorter cracks then filled in the 
pattern, as in Fig. |^a). For intermediate thicknesses 
of clay regular straight cracks appeared over each ridge, 
followed by cracks between these primary cracks, and at 
right angles to them, like the rungs of a ladder. As in 
Fig. gb), additional short cracks also often formed along 
the troughs of the pattern, especially in the later stages 
of drying. For thick layers, as in Fig. |^c), there was 
no obvious effect of the substrate on the crack pattern, 
which was essentially isotropic. 

Linear elastic fracture mechanics has no inherent 
length-scales, other than those set by the geometry of 
the system under study. If this assumption holds, then 
the dependence of the fracture pattern on the shape of 
our substrates can be entirely captured by considering 
two dimensionless parameters, such as the dimensionless 
layer thickness, h = iL/A, and a dimensionless ampli¬ 
tude, a = A/A, of the substrate’s relief. Since in many 
situations, including ours, the spacing of thin-film cracks 
scales with the thickness of the cracking layer 0 [13- 
lE], the parameter h should describe the relative match 
between the natural spacing of any cracks and the wave¬ 
length of the perturbations affecting their pattern. The 
dimensionless amplitude a describes the sharpness of the 
peaks on the substrate. 

For the linear wave plates we categorised each fully- 
dried crack pattern by visual inspection as either wavy, 
straight or isotropic. As summarised in Fig. the differ¬ 
ent patterns divide up the expected dimensionless phase 
space into well-defined regions. At the boundaries be¬ 
tween the different patterns were small areas of coex¬ 
istence, for example between wavy and straight cracks. 
This coexistence represents two scenarios - either both 
patterns were observed within different regions of a single 



FIG. 6: Example power spectra of ladder-like crack networks, 
with (a) cracks over each ridge of the substrate and (b) cracks 
over the ridges and troughs. Insets show representative sub¬ 
sections, 6 cm wide, of the original images that were summed, 
vertically, to give average crack densities along the x-axis, for 
Fourier processing. 


experiment, or were seen for the same h, but on different 
plates. For dimensionless layer heights in excess of about 

= 1, only isotropic patterns were ever seen. 

Although there are clear visual differences between the 
patterns shown in Fig. [^ there are no well-established 
methods by which they can be characterised. In order to 
analyse our patterns quantitatively we therefore need to 
develop new order parameters similar to, for example, the 
orientational order parameter of liquid crystals [37] . We 
present three such measurements next, based on Fourier 
methods, the relative orientation of the regions outlined 
by the crack network and the positions of neighbouring 
crack vertices. 


A. Fourier analysis of crack patterns 

The presence of evenly spaced cracks that are well- 
aligned with the substrate’s features suggests the use 
of Fourier methods to quantify our crack patterns. 
Although a two-dimensional transform could be used 
to measure crack spacings, we focus here on a one¬ 
dimensional approach, which allows us to measure the 
relative alignment and periodicity of the cracks with re¬ 
spect to the linear substrates. To this end, we first col¬ 
lapse each binary image to find the average density of 
cracks at each position along the x-direction. This is done 
by summing over the ^-direction of the images, where 
cracks are represented by ones and uncracked material 
by zeros. We then find the power spectrum of the crack 
density as the absolute value of the square of its Fourier 
transform, normalised to have an integrated power of one. 

In cases where cracks are aligned with the substrate 
there are sharp peaks in the power spectrum of the crack 
density when the wavenumber is an integer multiple 
of the wavelength A of the substrate, as in Fig. In 
other words, the substrate prefers there to be a inte¬ 
ger number of cracks, per wavelength. The responses 
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FIG. 7: Spectral power of the lowest-order peaks for cracks 
over plate 2. At /i = 0.6 there is a transition from a pattern 
with cracks along the peaks and troughs of the substrate, to 
one with cracks only along the peaks. Insets show 4x4 cm^ 
samples of crack patterns from the indicated conditions. 


of all five linear wave plates were very similar: typi¬ 
cally, most power was found in the first or the second 
peak, i.e. when = 1 or 2. Furthermore, there was 
a change in the dominant mode of the power spectra 
around h = 0.5 — 0.75. Below this, most power was in the 
= 2 peak, while above this the ^a,A = 1 peak was also 
strong. Visually, the point where these intensities cross 
over coincides with the transition between patterns with 
cracks along every ridge (Fig. I^a), ^a,A = 1), and pat¬ 
terns with cracks along both the ridges and troughs (Fig. 
[^b), ^a,A = 2). The appearance of wavy cracks, at very 
small /i, are associated with a decline in the power of the 
^a,A = 2 mode. Since they meander, their crack density 
is not as strongly localised as the ladder-like cracks, and 
the Fourier analysis reflects this. The isotropic cracks, 
around h = 1 or thicker, show no significant coupling to 
their substrate’s periodicity. 

An example of all of the above behaviour is given in 
Fig.Q for plate 2, where a = 0.25. The results for plate 
3 are consistent with this, and extend the range of sup¬ 
pressed ^a,A = 2 down to h = 0.2. For plates 4 and 5, 
where the substrate is relatively sharper (a = 0.5), the 
sequence of responses is identical, but shifted to slightly 
higher values of h: the transition from one to two cracks 
per wavelength is at h = 0.75, for example. For plate 
1, with shallow relief, this transition is shifted down to 
about 0.4-0.5, but the Fourier analysis does not produce 
very clear results, otherwise. These variations in the con¬ 
ditions leading to the different types of crack patterns are 
consistent with the visual classifications given in Fig. 

Thus, Fourier methods capture the regularity of the 
observed crack patterns, and quantify how the transi¬ 
tion from one to two cracks per wavelength occurs when 
the periodicity of the substrate is about 1.5 to 2 times 
the thickness of the cracking layer - i.e. approaching 
the natural crack spacing of the clay over a flat sub¬ 
strate. However, it does not provide a good indication 



FIG. 8: Ordering of cracked regions. Ellipses are fit to each 
connected region outlined by cracks, and the orientational 
ordering of the ellipses is then studied. Shown are the best- 
fit ellipses (red) with their major axes (green), overlaid on 
binary (white cracks on black background) representations 
of (a) well-ordered regions over a radial plate and (b) well- 
ordered and (c) isotropic regions over a linear plate. 


of the alignment of the various crack features with the 
substrate. The absolute intensity of the power spectral 
peaks, which could be such a measure, is too easily af¬ 
fected by noise to serve in this way. The next two sections 
aim to provide more robust measures of alignment. 


B. An orientational order parameter 


As the clay dries in our experiments, its cracks sub¬ 
divide and isolate individual regions of material, whose 
shapes are then determined by their sets of final bound¬ 
ing cracks. For example, in Fig. [^b) the cracks around 
the dried peds are mostly parallel or perpendicular to 
the substrate relief. The cracked regions therefore have 
rectangular shapes, and tile the plane like rows of bricks. 
In contrast, the cracks in Fig. [^c) are only weakly influ¬ 
enced by the substrate and define cracked regions with 
more arbitrary shapes and orientations. This section 
aims to measure the degree of alignment of these regions, 
with respect to the ripples along their substrate. 

From each image we extracted a binary skeleton of the 
crack network (see Sect. [hai . Each connected region 
in this image is identified and an ellipse is fit to cover 
its area. The major axis of this ellipse goes through the 
region’s centroid and its orientation, is chosen to min¬ 
imise the second moment of area of the region, around 
the major axis. Figure shows examples of elliptical fits 
of cracked regions with three different patterns. 

To characterise the crack-bounded regions, through 
their best-fit ellipses, we use an orientational order pa¬ 
rameter that measures the relative alignment of rod-like 
particles, or any objects that are invariant after rotation 
by 180°. This parameter is well-known in the physics of 
liquid crystals {e.g [37]), where it is used to describe the 
transition from an isotropic to a nematic (orientational 
order without positional order) phase, for example. In 
our case, for dimension d = 2, the general definition [38] 
of this nematic order parameter simplifies to 


= 


d{cos‘^ 0) — 1 
d — 1 


= (cos(26>)). 


( 1 ) 
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h=H/X 

FIG. 9: The orientational order parameter Si (Eq. de¬ 
scribes how well-aligned individual peds (continuous regions, 
entirely bounded by cracks) are with the orientation of their 
substrate. Results are given for wave plates with dimension¬ 
less amplitudes (a) a = 0.125, (b) a = 0.25 and (c) a = 0.5. 
The solid and dashed lines represent the mean and standard 
deviation of measurements of Si on flat plates. 


where the angled brackets represent averaging over all el¬ 
lipses. The orientation, of the major axis of each ellipse 
is measured with respect to the ^-axis for the linear wave 
plates, and the radial direction for the radially-symmetric 
substrates. 

The orientational ordering for all our experiments is 
shown in Fig. For each dimensionless amplitude a 
there is good agreement between the different experi¬ 
ments, although the radial plates give slightly higher val¬ 
ues of Si than the linear plates. For small enough /i, 
in particular around h = 0.5 and below, the peds are 
very strongly aligned with the relief of the substrate, and 
this alignment slowly decreases as the layer thickness in¬ 
creases, until Si becomes statistically equivalent to its 
value over the fiat plate above a relative thickness of at 
most h = 1.25. Furthermore, the intensity of the order¬ 
ing decreases with the relative sharpness of the plates, 
a, for experiments with the same h. However, unlike the 
Fourier analysis, or the order parameter we introduce in 
the next section to describe the relative alignment of the 
individual cracks, Si remains high for very low h. It is a 
sensitive measure of the straight-to-isotropic transition, 
but because it averages over the entire shape of a cracked 
region, it does not detect the meandering of the cracks. 


C. Crack alignment order parameter 

A second order parameter is introduced to describe 
the relative alignment of the cracks with their substrate. 
Working from the skeletons of the cracks, we first find 
the coordinates of all vertices, or intersections between 
cracks. These skeletons form closed networks, where all 
isolated features and spurs have been removed. The ver¬ 
tices are therefore points along the skeleton that touch at 
least three different cracked regions, and this property is 
used to identify them. Crack segments are then defined 
by the paths between any two neighbouring, connected 
vertices, and their orientation (j) by the vectors connecting 
these two endpoints. Since one end of a crack cannot be 
effectively distinguished from the other, the distribution 
of crack angles is collapsed onto the range from 0-180°. 
These angles are referenced to the y or radial directions, 
for the linear or radial wave plates, respectively. 

It should be emphasised that this method creates a 
very local definition of the orientation of a crack, and 
splits long cracks into many individual segments, between 
adjacent intersections with other cracks. This is different 
from some other network analyses of cracks [33l EH |40] , 
which focus instead on the entire length of a crack, re¬ 
gardless of any intersections it has with other features. 
An example of a crack pattern, and the relative abun¬ 
dance of crack orientations in it, is given in Fig. 

In our more well-ordered experiments, as in Fig. 
there are typically prominent cracks along the ridges of 
the substrate, and shorter perpendicular cracks that join 
these together. Their measure must reflect this symme¬ 
try. In analogy with Si, it should take a value of +1 for 
cracks that are either parallel or perpendicular to some 
reference direction, 0 for random patterns, and —1 for 
cracks that are all completely mz5-aligned with the ref¬ 
erence (here, cracks at ±45° to the substrate features). 
The simplest such definition is for a crack alignment or¬ 
der parameter 

S 2 = (cos(4</>)). (2) 

We calculated this parameter for each drying experiment, 
and summarise the results in Fig. 

We found that the crack alignment order parameter, 
£' 2 , matches the orientational order parameter of the 
cracked regions, Si, in most situations. There is also 
good data collapse between the linear and radial wave 
plates, and between all plates of the same relative ampli¬ 
tude a, confirming that our resuits do not depend on the 
absoiute scaie of the drying experiments. There is a ioss 
of ordering around h = 1, again at siightiy higher vaiues 
of h for the sharper substrate ridges, and siightiy iower 
vaiues for the shaiiower reiief. However, uniike there 
is a ciear downturn in S 2 at iow iayer thicknesses - the 
cracks become iess aiigned for very thin iayers. As it is 
a iocai measure, the crack orientation is more sensitive 
than our other metrics to any curvature of the cracks, 
and some of the decrease in S 2 for very iow h is due to 
the wavy cracks. Nevertheiess, in aii cases S 2 starts to 
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Crack angle 0 (°) 

FIG. 10: The distribution of crack orientations in a ladder¬ 
like pattern. The peak at 0° and 180° (these two angles being 
equivalent) reflects the dominant cracks parallel to the y-axis, 
while the peak at 90° comes from the shorter ‘rungs’ of the 
ladders. 



FIG. 11: The crack alignment order parameter S 2 (Eq. 
measures the relative alignment of the cracks with the sub¬ 
strate, and takes large values for networks with cracks that 
are parallel or perpendicular to the underlying relief. Shown 
are the results for (a) a = 0.125, (b) a = 0.25, and (c) a = 0.5. 
The solid and dashed lines represent the mean and standard 
deviation of aS '2 for the flat plate. 


fall well before we can see wavy cracks by eye. In Section 


IV we will show how this effect, and the other transitions 


in the crack pattern, can be explained by the relative 
forcing of the uneven substrate, which will be shown to 
be strongest at the intermediate layer thicknesses where 
both order parameters are largest. 


D. Extension to other crack networks 

The two order parameters Si and S 2 introduced above 
describe, respectively, the relative alignment of crack- 
bounded areas, or individual cracks, with respect to some 
known angle. In our experiments there was an obvious 
direction to choose for this angle, based on the substrates 
used - either the local radial vector, or one of the carte¬ 
sian axes. However these metrics can also be used to de¬ 
termine whether a crack network is aligned or not with¬ 
out knowing, in advance, any such preferred direction. 
By extending the liquid-crystal analogy of the order pa¬ 
rameter Si, one may use the local director, or average 
local orientation of the relevant objects, as the reference 
direction. In other words, one may first calculate the 
average 0 or <p, either over all cracks or areas in a re¬ 
gion of interest, or locally over neighbourhoods that are 
large enough to provide good statistical averaging. Then, 
in calculating the order parameters through Eqs. and 

all angles can be referenced to this local director ori¬ 
entation. This would allow one to compare the relative 
degree of alignment of natural crack patterns like polygo¬ 
nal terrain |4T1142] and the hydraulic channels of heavily 
fractured rock [39l |43], for example. It would also al¬ 
low one to better measure manufactured patterns such as 
craquelure [6l[15], or custom-designed micro-cracks [8l|9|. 


IV. MODELLING THE CRACK PATTERN 
TRANSITIONS 


When a crack opens, it releases stored strain energy 
from its environment. In Griffith’s model of fracture this 
energy is used to create the new surface area of the crack, 
as it grows [45] . In this section we consider a finite 
element model of a channel crack opening over a sinu¬ 
soidal surface, and use the Griffith criteria for fracture to 
predict where the first cracks will appear. We will show 
how the observed transitions from wavy to straight to 
isotropic cracks can be explained in terms of the relative 
strain energy released by a crack opening at different lo¬ 
cations in a thin elastic layer drying - or cooling, or with 
some other misfit strain - over such a rippled substrate. 

In particular, we consider a linear poroelastic [46] 
model for cracks in drying clay. The total stress a in 
the clay is a sum of the elastic stress a supported by the 
network of clay particles and the capillary pressure p of 
the interstitial liquid, such that 

^ij — ^ij P^ij • (3) 

The negative sign here accounts for the convention that a 
positive stress is tensile, while a positive pressure is com¬ 
pressive. The pressure is generated by the tiny curved 
menisci of the water still trapped in the clay, and for 
these conditions p < 0. The network stress can be given 
in terms of the strain e by the linear elastic constitutive 






































FIG. 12: Numerical solution of a channel crack opening over 
(a) a trough or (b) a peak of the substrate. The deformations 
caused by the cracks are exaggerated, so as to be visible. 
The greyscale shows how much of the pre-crack stress a*x is 
released by the crack - with black shown where all the stress is 
released, and white where the crack does not affect the stress 
at all. 


relationship 


E 

TT^ 



V 

1 - 2v 



( 4 ) 


where E is the Young’s modulus and v the Poisson’s ratio 
of the clay. The strain tensor is, in turn, derived from 
the elastic displacements, of the clay by 


= {piU^ + djUi)l2. ( 5 ) 

The pressure term in Eq. provides the driving force 
for fracture. As the clay dries p becomes more and more 
negative (its magnitude increases), requiring the clay to 
compensate with a tensile strain that can lead to crack 
growth. In order to instead model cooling, or tissue 
growth {e.g. for crocodile skin m or leaves [27]), this 
pressure can be replaced by an equivalent thermal con¬ 
traction term [23l [24] , or the effects of differential growth 
in a brittle layer and its substrate. 

We solved Eqs. [3]-[^ using Matlab’s finite-element rou¬ 
tines, on domains like those shown in Eig. namely 
semi-rectangular regions with straight edges at the top 
{z = 1, dimensionless units, scaled so that the aver¬ 
age film thickness H = 1), left {x = 0) and right 
{x = 10) sides, and a corrugated base where z = 1 — £ = 
Acos(27rx/A + -0). Here ip is an offset phase, allowing 
us to examine the strain energy released by fracture at 
various positions along the substrate, and £ is the local 
film thickness. Plane strain conditions, with no displace¬ 
ments allowed in the ^-direction, were used to model a 
channel crack advancing with a constant profile along 
the sinusoidal plate. The upper surface was taken to be 
traction free {axz = ^zz = 0 on 2 : = 1), while the lower 
surface was treated as no-slip, such that u = ^ there, 
and only vertical displacements were allowed {u^ = 0) 
along the right boundary. The left-face boundary was 
either treated like the right one, for an intact film, or 
as a traction-free boundary, for a cracked film. Einally, 
the pressure p was taken to be constant everywhere, as 
appropriate for quasi-static fracture. 


The work done by the growth of the channeling crack 
can be measured by integrating the product of the pre¬ 
crack stress, cr*, and the opening displacements of the 
crack, J, as in nulls]. This is, essentially, the sum of 
the forces on the incipient crack surfaces, acting over the 
distance they travel as the crack opens. The strain energy 
released, per unit area of new crack, is thus 

G = ^ ^ ^ SxKx + SzKzdz. (6) 


Eor a simulation with a specific A and A we calculated 
the pre-crack stress field once, using '0 = 0, and used 
this solution to find a* at different positions along the 
substrate. We then calculated the opening displacements 
of a crack at different '0. Accounting for the opening 
displacement to either side of a crack can be done by 
letting 5x = + Ux{-tp), and 6y = Uy{tp) - Uy{-'tp). 

This gives us the values of G for cracks growing along 
the substrate at different 0, ranging from the crest of 
the sine wave, to its trough. 

Eor a fiat substrate, this problem has been studied in 
depth by Beuth ED, Xia and Hutchinson [22] and Yin et 
al. mms]. They showed that G scales as 


~ 2 E ' 


( 7 ) 


where is a dimensionless parameter that only depends 
on the elastic mismatch between the cracking film and its 
substrate, and where E = E/{1 — and ao = —p{l — 
2z/)/(l — u). Eor example[2T], for z/ = 1/3, ^ = 0.71. We 
confirmed their solutions with respect to the exponential 
decay of the height-averaged stress away from the crack 
[ 22 ] and the dependence of g on u [21] and the depth of 
fracture Enim, and checked that g was independent of 
the average film thickness, the film’s elastic modulus E 
and the capillary pressure p. The mesh used, and the 
horizontal domain length, were found to be sufficient to 
agree with these known solutions to within 1%. 

Eor an uneven substrate G depends on location, mea¬ 
sured here by the offset 0. As the clay dries the relative 
shape of G{'ip) will not change, but rather its magnitude 
will simply increase along with the capillary pressure. 
The Griffith criteria predict that fractures can only grow 
once G reaches some critical value, Gc^ which is a ma¬ 
terial parameter (see e.g. [45]). The first cracks in the 
brittle drying layer should therefore form at the positions 
where the strain energy released by these cracks, z.e G, 
would be a maximum, as these are the first places where 
G = Gc- Eor a Poisson ratio of u = 1/3, we calculated 
how G varied with position in a wide range of experimen¬ 
tal conditions, from A = 0.01-0.90 and A = 0.5-8.0. To 
focus only on geometric effects we used Eq. to rescale 
G by the relative stiffness 2E/7r(jQ, so that G = £g. The 
results are summarised in Eig. [l^a), which shows the 
values of 0 that maximise G for various A and A. De¬ 
pending on the shape of the substrate and the thickness 
of the clay, the preferred positions of the earliest cracks 
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FIG. 13: Modelling the straight-to-wavy transition, (a) Pre¬ 
dicted positions of the first cracks, based on the offset where 
the energy release rate G is maximised. Red diamonds and 
blue circles show the experimental parameters where ladder¬ 
like or wavy cracks were seen, respectively (from Fig. |^. 
(b) Energy landscapes of (5(0) for three conditions, indicated 
by (1-3) in (a), which cross over the straight-to-wavy transi¬ 
tion. Under these conditions there is a bifurcation from one 
energy maximum for a crack on the top of a ridge, to two 
preferred locations, on either side of the ridge. 


FIG. 14: Modelling the straight-to-isotropic transition. For 
(a) the colour indicates the difference between the maximum 
and minimum values of G. The red circles have a radius 
proportional to the crack orientation order parameter, S 2 , for 
each experiment. The most well-ordered patterns correspond 
to the patch of high AG around the centre of the image. 
Panel (b) shows how S 2 is well-correlated with A(5, between 
experiments where straight or isotropic patterns were seen, 
and corresponding simulations. 


can be either (i) over the peaks of the substrate, (ii) over 
the troughs or (iii) at pairs of points displaced to either 
side of each peak. 

Figure a) also replots the observations from Fig. 
corresponding to the conditions where wavy and ladder¬ 
like cracks are found. For all experiments where straight 
cracks were seen the model correctly predicts that the 
first cracks should appear over the ridges of the substrate. 
Where wavy cracks are seen, the model instead gener¬ 
ally predicts two stable positions for cracks, along the 
sides of the ridges. Although our model has no dynamic 
component, this empirical agreement of the division of 
phase-space suggests that the wavy cracks result from 
cracks moving back and forth between their two bistable 
positions. The boundary between the straight and wavy 
cases corresponds to a forward pitchfork bifurcation, as 
sketched in Fig. (Hb). 

To consider the more gradual transition from ladder¬ 
like to isotropic cracks, we looked at the relative vari¬ 
ations in (5, defining AG as the difference between the 
maximum and minimum values of for each value of 

A and A studied. If AG is small, then there is only weak 
guidance for the cracks, and they should appear more 
random than in cases of larger variation of energy release 
rate. We tested this hypothesis by comparing AG with 


the crack alignment order parameter of our experiments. 
First, we calculated AG for each point in our simulation, 
mapping the results in Fig. [T^ alongside the data showing 
S 2 for our various experiments with linear wave plates. 
This figure shows that there is a neighbourhood of driving 
parameters, centred around X/H = 2.2 and A/H = 0.6, 
where AG is largest. Our experiments pass through this 
region, and show their most well-ordered patterns as they 
do so. To quantify this agreement, we then calculated 
AG for the conditions corresponding to each experiment 
where ladder-like or isotropic cracks were seen. As shown 
in Fig. mb), there is indeed a good correlation between 
the range of energy release rates allowed for any partic¬ 
ular geometry, and the relative alignment of the cracks 
that formed over it. 


Finally, we note that some insight into the limits of 
our results can be gained through a qualitative compar¬ 
ison to channel cracks in flat films and coatings. Here 
the stresses relaxed by the crack will decay away over a 
length comparable to the thickness of the film [20l |22] . 
This effect can also be seen in Fig. 12 For a wavy sub¬ 
strate this limit is approached when the amplitude of the 
relief goes to zero, and under such conditions isotropic 
cracks with a spacing ^ H are clearly expected. How¬ 
ever, for finite A/H the limit can still offer insight. For 
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very long wavelength ripples, large X/H^ each crack will 
see a local environment, of size ^ that is approxi¬ 
mately flat. Although these conditions were not probed 
by our experiments, we would also expect an isotropic 
pattern in this case, but one with a local modulation 
of the crack spacing reflecting the changes in local thick¬ 
ness. For the opposite limit, of very sharp features where 
X/H 1, the region around any crack will sample many 
ridges and troughs, and these local variations will tend 
to cancel each other out (see Saint-Venant’s principle e.g. 
@9]); again, isotropic cracks would be expected. In other 
words, for finite-amplitude relief a brittle layer must filter 
out the high and low-frequency information of its sub¬ 
strate, and any contraction cracks that form should re¬ 
flect this. 


V. CONCLUSIONS 

Our results explore the effects of an uneven sub¬ 
strate on the development of contraction cracks in thin 
films, and are inspired by natural patterns of cracks 
around buried craters, paintings, and vapour-deposited 
thin films. We showed that the wavy, periodic undulation 
of a substrate can cause the ordering of a crack pattern 
that forms above it, and has the strongest influence on 
this pattern when the thickness of the cracking layer is 
comparable to, or slightly smaller than, the periodicity 


of the substrate. 

Generally we found that for increasing layer thickness 
the crack pattern changed from one with wavy cracks 
running along the troughs of the substrate, to straight 
cracks appearing over each peak and trough of the sub¬ 
strate, to straight cracks on only the peaks of the sub¬ 
strate, and finally to a disordered or isotropic pattern. 
The wavy-to-straight transition was the most abrupt 
change, as evidenced by examples of single experiments 
that were wavy on one side of the plate, and straight on 
the other. All these patterns were quantified by Fourier 
methods and orientational order parameters, similar to 
the nematic order parameter of liquid crystals, which are 
easily adaptable to other types of crack networks. Fi¬ 
nally, we showed how the transitions between pattern 
type could be explained by a simple finite-element model 
of channeling cracks advancing over uneven relief, regard¬ 
less of whether those cracks were driven by drying, cool¬ 
ing or the differential growth of crocodilian skin. 
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